Unsupervised Learning with Linear Algebra

Reconstructing an image using Singular Value Decomposition (SVD)

Author

Anweshan Adhikari

Introduction

The use of unsupervised learning is not to make predictions, but to uncover hidden structures in the features of the data. In this blog post, we’re going to dive into one such unsupervised learning method known as Singular Value Decomposition (SVD).Singular Value Decomposition (SVD) is a fundamental technique in linear algebra specially in data compression and noise reduction. Also in this blog post we will explore an interesting application of SVD - image reconstruction. The SVD of a matrix \(A \in \mathbb{R}^{m \times n}\) is:

\(A=UDV^T\) ;

U and V are orthogonal matrices. This means that their columns are orthogonal unit vectors. D is a diagonal matrix containing singular values s\(\sigma_i\) of A. These values represent the magnitudeof each new dimension in the transformed data.

Converting image to Grayscale

Here I have chosen a RGB colored image of Mario and converted it to a grayscale image using the PIL package.

Code

from PIL import Image
from matplotlib import pyplot as plt 
import urllib
import numpy as np

def read_image(url):
    return np.array(Image.open(urllib.request.urlopen(url)))

url = "https://images-wixmp-ed30a86b8c4ca887773594c2.wixmp.com/f/6b51c512-9408-4f6c-9b9d-8bbddbafe45c/dfx6346-aa22997e-3fac-4663-88e8-581060c0dcf0.png/v1/fill/w_988,h_809/gcn_mario_kick__rgb__by_giltchyboi64_dfx6346-pre.png?token=eyJ0eXAiOiJKV1QiLCJhbGciOiJIUzI1NiJ9.eyJzdWIiOiJ1cm46YXBwOjdlMGQxODg5ODIyNjQzNzNhNWYwZDQxNWVhMGQyNmUwIiwiaXNzIjoidXJuOmFwcDo3ZTBkMTg4OTgyMjY0MzczYTVmMGQ0MTVlYTBkMjZlMCIsIm9iaiI6W1t7ImhlaWdodCI6Ijw9OTE5IiwicGF0aCI6IlwvZlwvNmI1MWM1MTItOTQwOC00ZjZjLTliOWQtOGJiZGRiYWZlNDVjXC9kZng2MzQ2LWFhMjI5OTdlLTNmYWMtNDY2My04OGU4LTU4MTA2MGMwZGNmMC5wbmciLCJ3aWR0aCI6Ijw9MTEyMyJ9XV0sImF1ZCI6WyJ1cm46c2VydmljZTppbWFnZS5vcGVyYXRpb25zIl19.Vb9e3dCybP27DsdJD3yf2MYC_M0MCm1B69gb9Wp0C-0"
img = read_image(url)

fig, axarr = plt.subplots(1, 2, figsize=(7, 3))

def to_greyscale(im):
    return 1 - np.dot(im[..., :3], [0.2989, 0.5870, 0.1140])

grey_img = to_greyscale(img)

axarr[0].imshow(img)
axarr[0].axis("off")
axarr[0].set(title="original")

axarr[1].imshow(grey_img, cmap="Greys")
axarr[1].axis("off")
axarr[1].set(title="greyscale")

[Text(0.5, 1.0, 'greyscale')]

Code

from PIL import Image
from matplotlib import pyplot as plt 
import urllib
import numpy as np

Let’s look at the shapes of the original image and grayscale image:

Code

a=np.shape(img)
b=np.shape (grey_img)
print(f"Shape of original image: {a}")
print(f"Shape of original image: {b}")

Shape of original image: (809, 988, 4)
Shape of original image: (809, 988)

The grayscale image retains the same width and height but loses the color channels. Now the grayscale image serves as a 2-dimensional matrix, which can be input into the SVD.

SVD Implementation

Code

def svd_reconstruct(image, k):
    # Performing SVD
    U, s, V = np.linalg.svd(image, full_matrices=True)
    D = np.zeros(image.shape)
    np.fill_diagonal(D, s)
    
    # Selecting
    D_k = D[:k, :k]
    U_k = U[:, :k]
    V_k = V[:k, :]

    # Reconstructing our image
    reconst_img = U_k @ D_k @ V_k
    
    return reconst_img

Code

k = 50 
reconst_img = svd_reconstruct(grey_img, k)

fig, axarr = plt.subplots(1, 2, figsize=(12, 6))

axarr[0].imshow(grey_img, cmap='Greys')
axarr[0].axis('off')
axarr[0].set_title('Original Image')

#reconstructing image using 50 singular values
axarr[1].imshow(reconst_img, cmap='Greys')
axarr[1].axis('off')
axarr[1].set_title('Reconstructed Image with k = ' + str(k))

plt.show()

At this k value, we can already see the image coming together which indicates that using only the 50 most important features can reconstruct our image to a decent level.

Experimentation

By setting k to 50, we are saying that we want to reconstruct our image using only the top 50 most important features. This compresses the image because fewer features are used to construct the image.

As k increases, the reconstructed image will become more and more similar to the original because more features are used in the reconstruction. Oppositely, As k decreases, less information is used in the reconstruction, and the image becomes more and more compressed. We can see this in the graph below

Code

import matplotlib.pyplot as plt

def storage_percentage(original, k):
    return (k * (original.shape[0] + original.shape[1]) + k) / (original.shape[0] * original.shape[1]) * 100

# Creating the subplots
fig, axarr = plt.subplots(2, 4, figsize=(20, 10))

#original iamge
axarr[0, 0].imshow(grey_img, cmap="Greys")
axarr[0, 0].axis("off")
axarr[0, 0].set_title("Original Image\n100% Storage")

#no. of features
k_values = [1,5, 25, 50,100, 175, 250]


for i, k in enumerate(k_values, start=1):
    # Reconstructing the image
    reconst_img = svd_reconstruct(grey_img, k)

    row = i // 4
    col = i % 4

    axarr[row, col].imshow(reconst_img, cmap="Greys")
    axarr[row, col].axis("off")
    axarr[row, col].set_title(f"Reconstructed Image\n k = {k}\n{storage_percentage(grey_img, k):.2f}% Storage")

plt.tight_layout()
plt.show()

As mentioned earlier, the selection of k is a balance between reducing storage space (lower k) and preserving image quality (higher k). It is interesting to see how efficiently an image can be compressed/reconstructed without losing much details from the original greyscale image. Moreover, there is a little observational difference when k is set to 175 and 250, however the percentatge of storage is significantly different. At the value of k =100, the image is easily recognizable while only taking 22.39% of the original storage space.

Optional Extras:

The compression factor is the percentage of singular values retained in the reconstruction. For example, a compression factor of 0.50 would mean retaining 50% of the original singular values. In the code below, we modified our svd_reconstruct method for users to specify a desired compression factor and select the number of components k to use based on this selection. I have demonstrated an example with a compression factor of 0.08

Code

def no_of_features(compression_factor, img):
    height, width = img.shape
    original_storage = height * width
    target_storage = compression_factor * original_storage
    k = target_storage / (height + width + 1)
    return round(k)

compression_factor = 0.08
reconst_img = svd_reconstruct(grey_img, no_of_features(compression_factor, grey_img))

Code

# Creating subplots
fig, axarr = plt.subplots(1, 2, figsize=(12, 6))

axarr[0].imshow(grey_img, cmap='Greys')
axarr[0].axis('off')
axarr[0].set(title='Original Image')

axarr[1].imshow(reconst_img, cmap='Greys')
axarr[1].axis('off')
axarr[1].set(title=f'Reconstructed Image\n k = {no_of_features(compression_factor, grey_img)}\n compression_factor= {compression_factor}')
plt.tight_layout()
plt.show()
#

Final Remarks

The use of SVD and the implementation of the svd_reconstruct function allow efficient image compression and uncover hidden structures in the features of the data. The ability to control the level of reconstruction could be a great tool for various image processing applications.